#clean up
rm(list=ls())

#packages
library(XLConnect)
library(irr)
library(tools)
library(gmodels)
library(maps)
library(foreign)
library(lattice)
library(plm)
library(lmtest)
library(car)
library(xtable)
library(systemfit)
library(stargazer)

#load data
#setwd('/Volumes/MONOGAN/immigration/stateLaws/')
#setwd('/Users/jamie/Documents/immigration/stateLaws/')
data<-read.csv("immLawsMasterClean.csv")
state.aggregate<-read.csv("stateAggregate.csv")
panel.data<-read.csv("panelImmig.csv")

###CREATE UNWEIGHTED MEASURES WITH AND WITHOUT RESOLUTIONS###
table(data$subject)
no.res<-subset(data,subset=subject!="resolution")

#Cross Section#
state.aggregate$pro.un<-as.matrix(by(as.numeric(data$pro==1),INDICES=data$state,FUN=sum))
#cbind(rownames(as.matrix(by(as.numeric(data$pro==1),INDICES=data$state,FUN=sum))),as.character(state.aggregate$state))
state.aggregate$con.un<-as.matrix(by(as.numeric(data$pro==0),INDICES=data$state,FUN=sum))
state.aggregate$score.un<-log((state.aggregate$pro.un+1)/(state.aggregate$con.un+1))

state.aggregate$pro.un.nores<-as.matrix(by(as.numeric(no.res$pro==1),INDICES=no.res$state,FUN=sum))
state.aggregate$con.un.nores<-as.matrix(by(as.numeric(no.res$pro==0),INDICES=no.res$state,FUN=sum))
state.aggregate$score.un.nores<-log((state.aggregate$pro.un.nores+1)/(state.aggregate$con.un.nores+1))

with(state.aggregate,cor(score,score.un))
with(state.aggregate,cor(score,score.un.nores))
with(state.aggregate,cor(score.un,score.un.nores))

#Panel Data#
data$state.year<-paste(data$state,data$year,sep=".")
no.res$state.year<-paste(no.res$state,no.res$year,sep=".")
panel.data$state.year<-paste(panel.data$state,panel.data$year,sep=".")

panel.pro.un<-as.data.frame(as.matrix(by(as.numeric(data$pro==1),INDICES=list(data$state.year),FUN=sum)))
panel.pro.un$state.year<-rownames(panel.pro.un)
colnames(panel.pro.un)[1]<-"pro.un"
panel.con.un<-as.data.frame(as.matrix(by(as.numeric(data$pro==0),INDICES=list(data$state.year),FUN=sum)))
panel.con.un$state.year<-rownames(panel.con.un)
colnames(panel.con.un)[1]<-"con.un"

panel.pro.un.nores<-as.data.frame(as.matrix(by(as.numeric(no.res$pro==1),INDICES=list(no.res$state.year),FUN=sum)))
panel.pro.un.nores$state.year<-rownames(panel.pro.un.nores)
colnames(panel.pro.un.nores)[1]<-"pro.un.nores"
panel.con.un.nores<-as.data.frame(as.matrix(by(as.numeric(no.res$pro==0),INDICES=list(no.res$state.year),FUN=sum)))
panel.con.un.nores$state.year<-rownames(panel.con.un.nores)
colnames(panel.con.un.nores)[1]<-"con.un.nores"

#merges
panel.data<-merge(x=panel.data,y=panel.pro.un,by="state.year",all.x=T)
panel.data<-merge(x=panel.data,y=panel.con.un,by="state.year",all.x=T)
panel.data<-merge(x=panel.data,y=panel.pro.un.nores,by="state.year",all.x=T)
panel.data<-merge(x=panel.data,y=panel.con.un.nores,by="state.year",all.x=T)

#deal with missingness
panel.data$con.un[is.na(panel.data$con.un)]<-0
panel.data$pro.un[is.na(panel.data$pro.un)]<-0
panel.data$con.un.nores[is.na(panel.data$con.un.nores)]<-0
panel.data$pro.un.nores[is.na(panel.data$pro.un.nores)]<-0
no.session<-panel.data$state.year%in%c('Montana.2006','Montana.2008','Montana.2010','Montana.2012','Montana.2014','Montana.2016','Nevada.2006','Nevada.2008','Nevada.2010','Nevada.2012','Nevada.2014','Nevada.2016','North Dakota.2006','North Dakota.2008','North Dakota.2010','North Dakota.2012','North Dakota.2014','North Dakota.2016','Texas.2006','Texas.2008','Texas.2010','Texas.2012','Texas.2014','Texas.2016','Oregon.2006','Oregon.2008','Oregon.2010','Arkansas.2006','Arkansas.2008')
panel.data$con.un[no.session]<-NA
panel.data$pro.un[no.session]<-NA
panel.data$con.un.nores[no.session]<-NA
panel.data$pro.un.nores[no.session]<-NA

#compute scores
panel.data$score.un<-log((panel.data$pro.un+1)/(panel.data$con.un+1))
panel.data$score.un.nores<-log((panel.data$pro.un.nores+1)/(panel.data$con.un.nores+1))

with(panel.data,cor(score,score.un,use="pairwise"))
with(panel.data,cor(score,score.un.nores,use="pairwise"))
with(panel.data,cor(score.un,score.un.nores,use="pairwise"))

###PANEL MODEL, UNWEIGHTED COUNTS###
panel.data$time<-panel.data$year-2005
panel.mod<-plm(score~cces.smooth+demUnif+repUnif+gsp+fb9015+fbPct2015+squireProfess+time,data=panel.data,model="random",index=c("state","year"));summary(panel.mod)
panel.mod.un<-plm(score.un~cces.smooth+demUnif+repUnif+gsp+fb9015+fbPct2015+squireProfess+time,data=panel.data,model="random",index=c("state","year"));summary(panel.mod.un)
panel.mod.un.nores<-plm(score.un.nores~cces.smooth+demUnif+repUnif+gsp+fb9015+fbPct2015+squireProfess+time,data=panel.data,model="random",index=c("state","year"));summary(panel.mod.un.nores)
stargazer(panel.mod,panel.mod.un,panel.mod.un.nores,covariate.labels=c("Public Opinion Ideology","Unified Democrat","Unified Republican","Per Capita GSP","Delta Foreign Born 1990-2015","Percent Foreign Born 2015","Legislative Professionalism","Time","Intercept"),star.cutoffs=c(.1),type="html")

###CROSS-SECTIONAL MODEL, UNWEIGHTED COUNTS###
cross.sec<-lm(score~cces.smooth+demUnif+repUnif+gsp+fb9015+fbPct2015+squireProfess,data=state.aggregate);summary(cross.sec)
cross.sec.un<-lm(score.un~cces.smooth+demUnif+repUnif+gsp+fb9015+fbPct2015+squireProfess,data=state.aggregate);summary(cross.sec.un)
cross.sec.un.nores<-lm(score.un.nores~cces.smooth+demUnif+repUnif+gsp+fb9015+fbPct2015+squireProfess,data=state.aggregate);summary(cross.sec.un.nores)
stargazer(cross.sec,cross.sec.un,cross.sec.un.nores,covariate.labels=c("Public Opinion Ideology","Unified Democrat","Unified Republican","Per Capita GSP","Delta Foreign Born 1990-2015","Percent Foreign Born 2015","Legislative Professionalism","Intercept"),star.cutoffs=c(.1),type="html")

